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Particle  size  distributions  in  atmospheric  clouds 

By  Roberto  Paoli  & Karim  Shariff  f 

In  this  note,  we  derive  a transport  equation  for  a spatially  integrated  distribution  function 
of  particles  size  that  is  suitable  for  sparse  particle  systems,  such  as  in  atmospheric  clouds. 
This  is  done  by  integrating  a Boltzmann  equation  for  a (local)  distribution  function  over 
an  arbitrary  but  finite  volume.  A methodology  for  evolving  the  moments  of  the  integrated 
distribution  is  presented.  These  moments  can  be  either  tracked  for  a finite  number  of 
discrete  populations  (“clusters”)  or  treated  as  continuum  variables. 


1.  Introduction 

Particles  are  present  in  atmospheric  clouds  in  several  forms  such  as  liquid  droplets,  non- 
volatile aerosols  or  ice  crystals.  Their  microphysical  properties  control  many  processes 
such  as  the  production  of  rain  in  stratocumulus  clouds  and  radiation  through  cirrus 
clouds.  These  properties  depend  on  the  way  particles  interact  with  the  surrounding  air, 
through  fluid-dynamic  and  thermodynamic  processes.  As  these  processes  usually  take 
place  at  small  spatial  scales,  the  interaction  of  particles  with  atmospheric  turbulence  is  an 
important,  though  complex,  problem  in  cloud  physics  (Shaw  2003).  Prom  a computational 
point  of  view,  two  major  factors  contribute  to  this  complexity.  First  is  the  very  high 
turbulence  Reynolds  number  and  the  large  range  of  spatial  scales  (Vaillancourt  & Yau 
2000;  Shaw  2003):  for  convective  clouds,  the  ratio  of  energy-containing  to  dissipative 
length  scales  is  O(105),  while  the  Reynolds  number  of  the  largest  eddies  is  0(1O6  to 
107).  The  second  factor  is  that  the  mean  distance  A between  particles  is  of  the  order  of 
the  Kolmogorov  scale  rj.  Thus,  if  one  contemplated  direct  numerical  simulation  (DNS) 
where  all  spatial  scales  are  resolved,  then  one  would  have  to  track  individual  particles. 
Since  DNS  resolution  is  not  affordable  for  these  flows,  Eulerian  formulations  for  the 
liquid/solid  phase  are  widely  used  in  the  simulation  of  clouds.  These  formulations  fall 
into  two  main  classes.  The  first  is  a “two  fluid  model”  where  particles  are  modeled  as 
a continuum  having  a local  mass  density  per  unit  volume.  This  approach  carries  no 
information  about  the  distribution  of  the  particle  size.  In  the  second  approach,  some 
physical  properties  (e.g.  the  mean  radius)  of  some  “ensemble”  of  particles  are  explicitly 
computed  at  each  physical  location  x.  The  concept  of  particle  size  distribution  at  a 
point  at  this  stage  of  our  discussion  is  ambiguous  but  will  be  clarified  later.  A standard 
procedure  used  in  two-phase  flow  models  (e.g.  Williams  1965;  Cotton  & Anthes  1989; 
Crowe  et  al.  1998)  to  describe  an  ensemble  of  particles  is  to  define  a distribution  function 
/,  in  a manner  analogous  to  the  kinetic  theory  of  gases.  In  kinetic  theory,  a distribution 
function  /(x,  v;t)  is  defined  where  /(x,  v;  t)  £x  <5v  represents  the  number  of  molecules 
that  at  time  t are  between  x and  x + 5x  and  whose  velocity  is  between  v and  v + Sv. 
It  is  assumed  that,  as  (<5x,  <Jv)  — + 0,  the  phase  volume  still  contains  a sufficently  large 
population  of  molecules  that  statistics  can  be  used.  This  is  usually  true  in  gasdynamics 
because  the  mean  free  path  of  molecules  is  much  smaller  than  the  continuum  scale  one 
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cares  about.  The  extension  of  this  appraoch  to  particles  other  than  molecules  is  formally 
straightforward  (see  for  example  the  book  by  Williams  1965),  as  long  as  the  continuum 
description  remains  valid.  In  the  case  of  atmospheric  clouds  however,  A ~ 77  so  that  only 
a few  particles  rather  than  a population  are  present  in  a volume  V = 0(j ?3).  The  object 
of  this  note  is  first  to  derive  a transport  equation  for  an  integrated  distribution  function 
To,  describing  an  ensemble  of  particles  inside  an  arbitrary  but  finite  volume  Vo-  Then, 
this  approach  is  specialized  to  atmopsheric  clouds.  Finally,  a methodology  is  proposed 
to  solve  for  the  moments  of  the  distribution  function  To- 


2.  Distribution  function  in  atmospheric  clouds 

Let  us  consider  a population  of  particles  in  a cloud  from  an  ensemble  of  realizations. 
At  any  time  t,  each  particle  p occupies  the  position  xp(t)  in  physical  space,  moves  with 
velocity  u p(t)  and  changes  its  radius  rp(t).  This  population  can  be  represented  in  a 
phase  space,  defined  by  the  generalized  coordinates  q(f)  and  evolving  via  the  generalized 
velocities  U(q(f)): 


q(t)  = 


' r(t)  ' 

x(t)  , 

. UW  . 


u«>«»  = § 


r(q(t)) 

u(q(t)) 

F(q(0) 


(2.1) 


where  x and  u are  the  spatial  coordinates  and  velocities,  F is  the  functional  law  of  the 
force  per  unit  mass  acting  on  the  particle  and  r is  the  functional  law  of  the  growth  rate 
of  its  radius.  At  any  time  t,  each  particle  of  the  population  occupies  a point  q(f)  in  this 
space.  Instead  of  tracking  each  particle  we  wish  to  follow  the  evolution  of  a distribution 
function  /(q(f);£)  of  the  population.  This  is  defined  in  such  a way  that 


/(q(f);f)  5Q{t) 


(2.2) 


is  the  number  of  particles  that  at  time  t are  inside  a cube  of  volume  6Q(t)  in  phase  space, 
located  between  the  coordinates  q(f)  and  q(f)  + Sq(t).  After  a time  dt,  this  volume  has  a 
value  SQ(t  + dt)  and  the  diagonally  opposite  corners  of  the  cube  are  mapped  to  q(t  + dt), 
and  q(t  + dt)  + 5q(f  + dt),  respectively.  At  the  same  time,  each  particle  can  change  its 
q and  the  particle  number  can  vary  because  of  evaporation  or  coagulation.  To  derive 
a transport  equation  for  / we  need  to  relate  these  quantities.  For  the  sake  of  clarity, 
we  will  drop  the  explicit  dependence  on  t in  all  variables,  and  define  t'  = t + dt,  and 
q'  = q(f  + dt).  Then  we  have 

q7  = q + U(q)  dt  (2.3) 

q'  + <5q'  = q + <5q  + U(q  + Sq)dt  (2.4) 


The  last  term  in  (2.4)  can  be  expanded  as  U(q  + <Sq)  = U(q)  + (VqU)<5q  where  VqU  is 
the  gradient  of  U in  phase  space.  Substituting  (2.3)  into  (2.4)  one  gets 

<5q'  = (I  + dfVqU)<5q  (2.5) 

The  change  of  phase  space  volume  5Q'  —SQ  is  related  to  the  divergence  of  the  generalized 
velocity  U by 


SQ'  - SQ 


.2 
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Let  K (q;  t)  describe  the  general  rate  of  gain  or  loss  of  particle  number  due  to  coagulation 
or  evaporation, 

/( q';  t')  SQ'  = /( q;  t)  6Q  + K{ q;  t)dt  SQ  (2.7) 

Substituting  (2.3),  (2.4)  and  (2.6)  into  (2.7)  and  expanding  gives: 


/+(§£  + V^u)*  + 0<^ 


1 + Vq-U  dt  + 0{dt2) 


SQ  = fSQ  + KdtSQ  (2.8) 


where  all  quantitites  are  evaluated  at  t.  Taking  the  limit  dt  — > 0 and  neglecting  higher 
order  infinitesimals  gives  a Boltzmann  equation 


^ + Vq.(/U)  = tf  (2.9) 

Finally,  inserting  the  different  components  of  q and  U by  means  of  (2.1)  gives  the  more 
usual  form  employed  in  two-phase  flow  literature  (e.g.  Williams  1965): 


^ + Vx  • (/  u)  + Vu  • (/  F)  + = K (2.10) 

The  term  F in  (2.10)  is  the  aerodynamic  drag  induced  by  the  flow  on  particles,  F « 
(up  — u /)/tp  (see  Crowe  et  al.  1998)  where  u / is  the  fluid  velocity  and  rp  = Apprp/18fj, 
is  a relaxation  time.  If  the  size  of  the  particle  rp  is  small,  rp  is  also  small  and  the  particle 
velocity  immediately  adjusts  to  the  flow  velocity.  In  the  following,  we  restrict  our  analysis 
to  this  case,  so  there  is  no  dependence  on  F in  (2.10): 


d£ 

dt 


+ Vx  • (/  u)  + 


9(fr) 

dr 


= K 


(2.11) 


3.  Integrated  distribution 

We  now  derive  an  integrated  version  of  (2.11).  Consider  a point  x0  in  physical  space  and 
an  arbitrary  (but  finite)  volume  Vo(xo(f)  around  it.  Then,  a space  integrated  distribution 
function  .Fo(r,Xo;f)  can  be  defined  as 

•F0(r,x0;t)=  f f(r,x;t)dV0(x)  (3.1) 

JVo 

so  that  3~o{r,  xo;  t)  Sr  represents  the  number  of  particles  that  at  time  t are  inside  a finite 
volume  V0  around  xo  and  whose  radius  is  in  between  r and  r + Sr  (the  dimensions  of  this 
function  are  [To[  = L~l  whereas  [/]  = L~l  x L~3).  The  mean  ipo  within  the  volume  Vo 
of  any  quantity  ij}(r,  x;  t)  associated  with  each  particle  is 


■0o(r,xo;f)  = — — ^ — [ f(r,x;t)ip(r,x-,t)dV0(x)  (3.2) 

J-o(r,Xo)t)  Jy0 

Let  us  introduce  a local  coordinate  y around  Xo,  x = x(x0,  y)  = Xo  + y,  so  that 

Vx-(.)  = VX0-(.)  + Vy-(.).  (3.3) 

Using  (3.2)  for  u and  r,  their  integral  values  over  Vq  become 


u0(r,xo;t)  = 


7-7  / f{r,  x0  + y;  t)  u(r,  x0  + y;  t)  dV0(y) 
To {r,x0;t)  JVo 

77~— A / f(r,x0  + y;t)r(r,x0  + y;t)dV0(y) 
j-o(r,x0,t)  jVo 


(3.4) 


r0(r,x0',t)  = 


(3.5) 
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Using  (3.3)  to  express  the  divergence  term  Vx  • (/  u)  = VX0  ■(/u)  + Vy'(/u),  a transport 
equation  for  To  can  be  derived  by  integrating  (2.11)  inside  volume  Vo: 

jv  %dV0  + VX0  • (£  fudVo^j  + £ Vy  • (/u)  dVo  + 

KdV0  (3.6) 

As  in  general  Xo  and  Vo  vary  in  time,  we  need  to  switch  volume  integration  and  time 
derivative  in  the  first  term  of  the  above  equation.  Using  Leibnitz  rule  and  (3.1),  one  has: 

TT  s I U dv°  - i % « + i ds°  W 

where  So  is  the  surface  around  Vo  and  S/  is  the  velocity  of  Sq  with  respect  to  the  fixed 
reference  frame.  Using  (3.7),  Gauss  theorem  and  the  definitions  (3.4)  and  (3.5)  in  (3.6), 
and  introducing  the  relative  velocity  with  respect  to  S0,  w = u — S/,  finally  gives 

^ + VX0  • (fo«o)  + + £ / w • n0  dS0  = K0  (3.8) 

which  is  a Boltzmann  equation  for  the  integral  distribution  function  To-  Equation  (3.8) 
formally  differs  from  (2.11)  because  of  the  surface  integral  in  the  left-hand  side.  This 
contains  the  (unknown)  local  distribution  function  /,  which  must  be  modeled  in  some 
way.  Note,  however,  that  if  Vo  is  a material  volume,  then  u = S/  on  5b,  and  the  surface 
flux  goes  to  zero. 


Ergodic  hypothesis 

In  the  previous  derivation  we  had  to  introduce  an  ensemble  of  realizations  in  order  to 
derive  a local  Boltzmann  equation  which  we  then  integrated  over  a finite  volume.  We 
now  make  the  hypothesis  that  Vo  is  large  enough  to  contain  a population  of  particles 
that  we  can  by-pass  the  ensemble.  In  other  words,  we  hypothesize  that  the  integrated 
Boltzmann  equation  (3.8)  is  valid  for  a single  realization  if  Vo  is  large  enough.  Typically, 
the  grid  size  in  cloud  codes  is  A > 1 m while  A ~ tj  « 10~3  m,  so  each  grid  cell  contains 
at  least  109  particles. 


3.1.  Ensemble  averages 

The  total  number  of  particles  No  in  spatial  volume  Vo  can  be  obtained  by  integrating  To 
over  all  possible  radii, 

POO 

JV0(x0;t)  = / T0{r,x0;t)dr.  (3.9) 

Jo 


Using  (3.9),  the  ensemble  average  (^o)  of  any  variable  ip  is  obtained  by  integrating  (3.2) 
over  r: 


1 f°° 

{ip0)(x0;t)  = No(Xo.tj  J ipo(r,Xo;t)To(r,Xo]t)dr 

In  particular,  the  ensemble  velocity  and  radius  growth  rate  are 

1 f°° 

(u0)  (x0;  t)  = ^ u0(r,  x0;  t)  T0{r , x0;  t)  dr, 

1 f°° 

(r0)  (xq;  t)  = N(./xLyt<;  Jo  ro(r,  x0;  t)  To{r,  x0;  t)  dr 


(3.10) 


(3.11) 
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For  further  analysis,  it  will  be  useful  to  introduce  the  mean  radius,  (ro),  and  the  variance, 
(A rg),  of  the  population.  Using  (3.2),  these  are  given  by 

/OO 

rT0(r,x0;t)dr  (3.12) 

rOO 

(A rg)(x0;f)  = lr  ~ (ro)  (u0;  t)f  T0{r,x0\t)  dr  (3.13) 

Note  that  u(r,  y;t)  and  (u0)  (xo;  t)  (same  for  r and  (ro))  have  a different  physical  mean- 
ing: u(y)  represents  the  velocity  of  a particle  in  a neighboorood  of  y.  On  the  other  hand, 
Uo(xo)  represents  a statistical  average  within  a population  of  particles  and  is  a continuum 
velocity  field,  associated  to  any  point,  xo,  of  the  physical  domain.  Indeed,  one  could,  in 
principle,  obtain  (uq)  and  (ro)  as 


1 Wo  i No 

<u°>  = E u*>  > 


P=  1 


No 


(3.14) 


p=i 


where  up  and  rp  are  the  velocity  and  the  radius  growth  rate  of  particle  p inside  V0.  In 
the  limit  of  No  — ■»  oo,  (3.11)  and  (3.14)  are  equivalent  but  we  only  have  access  to  To 
because  the  details  concerning  up  and  rp  of  each  physical  particle  are  unknown. 


The  next  step  is  to  relate  the  continuum  fields  (uo)  and  (r0)  to  the  corresponding  flow 
variables.  As  we  do  not  consider  here  any  force  acting  on  particles  (Sec.  2),  they  are 
simply  convected  by  the  fluid.  Therefore,  there  is  no  reason  why  two  particles  of  the 
same  population  and  different  radius  should  have  different  velocities,  i.e.  u0  is  statistically 
uncorrelated  with  r, 

u 0(r,x0;f)  = (u0)  (x0;t)  = u/(x0;t)  (3.15) 

where  u/(xo;  t)  is  the  fluid  velocity  at  x0.  The  same  arguments,  cannot  be  applied  to  ro, 
i.e.  ro  7^  (r)0,  because  each  particle  of  the  population  may  have  a different  growth  rate 
due  to  different  ’’reactions”  to  turbulent  fluctuations  in  the  flow-field,  as  discussed  next. 


4.  Particle  growth  by  condensation 

The  growth  of  the  radius  of  a single  particle  in  a medium  at  rest  can  be  simply  derived 
by  considering  a diffusion  equation  for  water  vapor  on  a particle  surface  (Pruppacher  & 
Klett  (1997)  p.  502)  and  is  given  by 

dr  D (pv  - psv(T))  PS 

dt  Tpw  r Tpw  r ^ ' 

where  D is  the  diffusion  coefficient  of  water  vapor  in  air  and  T is  the  psychrometric 
correction  associated  with  the  latent  heat  of  condensation;  and  pw  is  the  density  of  ei- 
ther water  or  ice.  The  vapor  densities  pv  and  psv{T)  are  evaluated,  respectively,  at  some 
“ambient”  condition  far  from  the  particle  and  at  the  surface  of  the  particle  (which  co- 
incides with  the  saturation  value  because  vapor  there  is  in  thermodynamic  equilibrium 
with  water/ice).  Thus,  in  such  a single-particle  picture,  the  radius  growth  rate  is  only 
controlled  by  the  supersaturation  S = pv  — psv.  As  first  pointed  out  by  Srivastava  (1989) 
(see  also  Khvorostyanov  & Curry  1999),  this  description  cannot  be  extended  straightfor- 
wardly to  a population  of  particles.  In  fact,  even  in  the  absence  of  turbulence  and  uniform 
S initially,  the  available  vapor  in  a cloud  is  not  equally  distributed  among  all  particles 
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because  of  their  random  spatial  distribution,  so  that  the  effective  supersaturation  avail- 
able at  a droplet  surface  (“microsaturation”),  can  differ  significantly  from  the  overall 
ensemble  averaged  supersaturation  (“macrosaturation”)  (Srivastava  1989).  In  addition, 
in  a turbulent  cloud,  each  particle  “reacts”  in  a different  way  to  turbulent  fluctuations  in 
the  flow-field:  for  example,  if  a supersaturation  fluctuation  arises,  it  will  be  absorbed  by 
each  particle  through  a complex  diffusional  process  of  vapor  involving  all  elements  of  the 
population  (Srivastava  1989).  Several  approaches  have  been  developed  in  the  atmospheric 
science  literature  (Pruppacher  & Klett  1997)  to  try  to  solve  this  complex  problem.  One 
of  these,  the  so-called  “stochastic  condensation”  approach  performs  Reynolds  averaging 
on  the  equation  for  condensational  growth,  resulting  in  covariances  that  can  be  though 
as  “Reynolds  stresses”  (Shaw  2003).  In  particular,  we  follow  Khvorostyanov  & Curry 
(1999)  (see  also  Pruppacher  & Klett  (1997)  p.505)  who  use  kinetic  theory  to  relate  the 
micro-  and  macro-saturation  in  a cloud.  Their  arguments  are  as  follows.  For  the  moment 
consider  the  situation  where  the  particle  radius  is  so  small  to  be  comparable  with  the 
mean  free  path  of  vapor  molecules.  In  this  case,  one  should  account  for  the  Brownian 
motion  of  molecules,  that  is  the  diffusion  associated  to  the  (random)  molecular  impact  on 
particles  surface.  As  shown  by  Pruppacher  & Klett  (1997),  this  can  be  done  in  (4.1)  by 
introducing  a modified  diffusion  coefficient  D*(r)  which  depends  linearly  on  the  radius  r 
(Pruppacher  & Klett  (1997)  p.506).  It  can  be  argued  (Crowe  et  al.  1998;  Khvorostyanov 
& Curry  1999)  that  the  effects  of  turbulent  fluctuations  of  vapor  density  or  in  super- 
saturation is  similar  to  Brownian  motion,  whereas  the  molecular  impact  on  particles  is 
substituted  by  their  interaction  with  turbulent  eddies  (note  that  this  picture  can  also 
be  extended  to  account  for  equivalent  Brownian  dispersion  of  particles,  induced  by  fluc- 
tuating fluid  forces  rather  than  density  fluctuations)  (Crowe  et  al.  1998).  The  “micro” 
supersaturation  Sp  available  to  particle  p of  an  ensemble  is 


Sp  = 


(4.2) 


where  (So)  (xo;t)  is  the  ensemble  supersaturation  available  to  the  population  within 
volume  Vq.  It  represents  the  supersaturation  that  would  be  at  Xo  if  there  were  no  particles, 
then  it  can  be  thught  as  the  fluid  supersaturation  at  Xo,  (So)  (x0;  t)  = 5/(xo;t).  Using 
(4.2)  and  the  previous  formalism  ( rp  — > r;  Sp  — ► 5(r,  y);  rp  —>  f (r,  y)),  one  gets  to 


S(r,  y)  _ (Sp)  _ Sf 
r (r0)  (ro) 

r{r,y ) = ggfor>=.  DSJL 
Vpwr  Tpw  (r0) 


Substituting  (4.3)  and  (4.4)  into  (3.11)  finally  gives 


(ro) 


1 

Aq(x0;  t) 


J-q  fo  dr 


DS(r,  y) 
r Pw r 


dV0  dr  = 


D Sf 
P Pw  (ro) 


(4.3) 

(4.4) 


(4.5) 


5.  Method  of  moments 

Even  neglecting  the  surface  term,  (3.8)  is  a p.d.e.  in  four-dimensional  space  (r,Xo;t) 
that  can  only  be  solved  numerically.  Some  atmospheric  cloud  codes  solve  a transport 
equation  for  a distribution  function  by  discretizing  the  particle  size  r in  a finite  number  of 
bins,  at  each  grid  location  (although  it  is  not  explicitly  mentioned,  they  are  conceptually 
discretizing  (3.8)  with  neglected  surface  terms). 
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In  this  section,  we  present  a simulation  strategy  based  on  the  method  of  moments 
proposed  in  Paoli  et  al.  (2002).  The  moments  mk  of  the  distribution  To  are  defined  by 


mk(xQ;t ) - 

Jo 


rkT0(r,x0;t)dr 


(5.1) 


Multiplying  (3.8)  by  rk  gives 


s or  rfcjr° dr ) + Vx°  ■ (/°°  rkjr°  u° dr) = +k  si  rk~ijr^° dr  (5-2) 

Using  (3.15),  (4.5)  and  (5.1)  and  assuming  that  To  —■ > 0 sufficiently  fast  as  r — > oc,  (5.2) 
becomes 

^i  + V„-(o  = *">*-'  (5.3) 

Under  all  assumptions  made,  (5.3)  describes  the  evolution  of  the  moments  of  the  integral 
distribution  function  To-  An  attractive  property,  deriving  from  the  microsaturation  model 
(4.2),  is  that  the  evolution  of  the  kth  moment  only  depends  on  the  previous  order  moment 
which  allows  one  to  close  the  system  (5.3)  without  any  further  assumptions  and  without 
presuming  the  shape  of  To-  Using  (3.9)-(3.13)  and  (5.1),  the  zero  and  the  first  two 
moments  are  easily  found  and  are  related  to  the  ensemble  average  radius  and  variance, 


too  = No  , mi  = N0  ( r0 ) , 


m2  = N0 


<A rS)  + (rofl  . 


(5.4) 


The  corresponding  evolution  equations  are  (we  put  a = D/Tpw) 


^ + VXo-(U/lVo)  = 0 
dm-  + VX0  • (uf  mi)  = aSf 


dt 

dm2 

dt 


+ VXo  • (Uf  m2)  = 2aSfN0 


(5.5) 

(5.6) 

(5.7) 


These  equation  are  coupled  to  the  continuum  fluid  phase  through  u/  and  Sf.  In  partic- 
ular, an  increase  in  particle  radius  by  condensation  implies  vapor  depletion  pv,  with 


• • [ a 2 • r j 4 1TpwD  Sf  f 2 t j a a 

Pv  = - I ATrpwr^r0 To  dr  = T / r‘T0dr  — -4npwaSf 

JO  1 Pw  VO)  Jo 


Nom2 

mi 


(5.8) 


The  usual  convection-diffusion-reaction  equation  for  a scalar  Yv  in  conservative  form 
(pv  = Yvpf  where  pf  is  the  total  gas  phase  density)  then  becomes: 

^ + VX0  • (uf  Pv)  + VX0  • (pfD  VXoyv)  = pv  = -4tt PwaSf  ^ (5.9) 

ox  m i 


Under  all  approximations  made,  (5.5)-(5.7)  and  (5.9)  (together  with  Navier-Stokes  equa- 
tions) describe  the  evolution  of  the  first  moments  of  the  size  distribution  of  a population 
of  particles.  These  moments  can  be  solved  by  using  either  an  Eulerian  or  a Lagrangian 
description  as  discussed  next. 


Eulerian  description 

In  this  case,  one  has  to  solve  for  mk  using  (5.5)-(5.7)  with  the  further  condition  that  the 
volume  V0  is  constant  in  time  (in  discretized  form  it  can  be  the  volume  of  a grid  cell).  It  is 
worth  mentioning  that  To  and  mk  are  continuous  functions  of  space  (not  grid  averages!), 
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so  their  gradients  contain  all  spatial  fluctuations  in  a turbulent  flow.  In  particular,  if  they 
are  filtered  in  a LES  approach,  the  correlations  rnpaj  and  m*,  5/  exist  at  subgrid  scale 
level  and  must  be  modeled.  ' 


An  approximate  Lagrangian  method 

Let  us  devide  the  total  number  of  particles  in  the  cloud  into  Nc  ” clusters” , each  containing 
a fixed  number,  Nj,  of  particles  (where  j = 1 ,..,NC).  The  position  of  each  cluster  is 
assumed  to  evolve  according  to 

^ = U ,(4)  (5-10) 

where  xj  is  the  center  of  the  volume  V0j  containing  cluster  j and  u/(xj)  is  the  fluid 
velocity  at  Xq.  Note  that  we  are  assuming  that  the  cluster  advects  rigidly  without  de- 
forming. Introducing  the  total  derivative  d()/dt  — dQ/dt  + VXo()u/  in  (5.5)-(5.7),  one 
can  “track”  the  moments  of  each  cluster  j as  (note  that  the  zeroth  moment  equation, 
Nj  = const  is  now  trivial) 


dm\ 

dt 


= aSf- 


N? 


dt 


= 2aSfNj 


(5.11) 

(5.12) 


where  Sf  = S/(x£)  is  the  fluid  supersaturation  at  xj.  The  advantage  of  the  Lagrangian 
approach  is  that  the  surface  term  in  (3.8),  which  was  neglected  to  get  to  (5.5)-(5.7),  is 
now  zero  because  V$  is  a material  volume,  u = on  S%  for  all  j. 


6.  Conclusions 

In  this  note  we  derived  a transport  equation  for  the  radius  distribution  function  of  a 
population  of  particles  in  an  atmospheric  cloud.  We  used  a simple  stochastic  condensation 
model  for  the  radius  growth  (taken  from  the  atmospheric  science  literature)  to  relate  the 
microsaturation  around  each  particle  to  the  macrosaturation  of  the  entire  population. 
Finally,  we  described  a procedure  to  solve  for  the  moments  of  the  distribution,  and  showed 
that  this  can  be  either  used  in  Eulerian  continuum  formulations  or  Lagrangian  tracking 
of  “clusters”  of  particles.  Future  plans  include  testing  of  this  method  by  comparison 
with  DNS  of  homogeneus  and  isotropic  turbulence  and  individual  particle  tracking;  and 
application  of  the  method  to  natural  and  contrail-generated  cirrus. 
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